home *** CD-ROM | disk | FTP | other *** search
/ Language/OS - Multiplatform Resource Library / LANGUAGE OS.iso / cpp_libs / newmat03.lha / newmat03 / newmat5.cxx < prev    next >
C/C++ Source or Header  |  1993-08-08  |  8KB  |  283 lines
  1. //$$ newmat5.cxx         Transpose, evaluate etc
  2.  
  3. // Copyright (C) 1991: R B Davies and DSIR
  4.  
  5. #include "include.hxx"
  6.  
  7. #include "newmat.hxx"
  8. #include "newmatrc.hxx"
  9.  
  10. //#define REPORT { static ExeCounter ExeCount(__LINE__,5); ExeCount++; }
  11.  
  12. #define REPORT {}
  13.  
  14.  
  15. /************************ carry out operations ******************************/
  16.  
  17.  
  18. GeneralMatrix* GeneralMatrix::Transpose(MatrixType mt)
  19. {
  20.    if ((int)mt==0) mt = Type().t();           // type of transposed matrix
  21.    GeneralMatrix* gm1 = mt.New(ncols,nrows); gm1->ReleaseAndDelete();
  22.  
  23.    if (mt == Type().t())
  24.    {
  25.       REPORT
  26.       for (int i=0; i<ncols; i++)
  27.       {
  28.      MatrixRow mr(gm1, StoreOnExit+DirectPart, i);
  29.          MatrixCol mc(this, mr.Store(), LoadOnEntry, i);
  30.       }
  31.    }
  32.    else
  33.    {
  34.       REPORT
  35.       MatrixRow mr(this, LoadOnEntry);
  36.       MatrixCol mc(gm1, StoreOnExit+DirectPart);
  37.       int i = nrows;
  38.       while (i--) { mc.Copy(mr); mr.Next(); mc.Next(); }
  39.    }
  40.    tDelete(); return gm1;
  41. }
  42.  
  43. GeneralMatrix* SymmetricMatrix::Transpose(MatrixType mt)
  44. { REPORT  return Evaluate(mt); }
  45.  
  46.  
  47. GeneralMatrix* DiagonalMatrix::Transpose(MatrixType mt)
  48. { REPORT return Evaluate(mt); }
  49.  
  50. BOOL GeneralMatrix::IsZero() const
  51. {
  52.    REPORT
  53.    real* s=store; int i=storage;
  54.    while (i--) { if (*s++) return FALSE; }
  55.    return TRUE;
  56. }
  57.  
  58. GeneralMatrix* ColumnVector::Transpose(MatrixType mt)
  59. {
  60.    REPORT
  61.    GeneralMatrix* gmx = new RowVector; MatrixErrorNoSpace(gmx);
  62.    gmx->nrows = 1; gmx->ncols = gmx->storage = storage;
  63.    return BorrowStore(gmx,mt);
  64. }
  65.  
  66. GeneralMatrix* RowVector::Transpose(MatrixType mt)
  67. {
  68.    REPORT
  69.    GeneralMatrix* gmx = new ColumnVector; MatrixErrorNoSpace(gmx);
  70.    gmx->ncols = 1; gmx->nrows = gmx->storage = storage;
  71.    return BorrowStore(gmx,mt);
  72. }
  73.  
  74. GeneralMatrix* GeneralMatrix::Evaluate(MatrixType mt)
  75. {
  76.    if ((int)mt==0) { REPORT return this; }
  77.    if (mt==this->Type()) { REPORT return this; }
  78.    REPORT
  79.    GeneralMatrix* gmx = mt.New(nrows,ncols);
  80.    MatrixRow mr(this, LoadOnEntry);
  81.    MatrixRow mrx(gmx, StoreOnExit+DirectPart);
  82.    int i=nrows;
  83.    while (i--) { mrx.Copy(mr); mrx.Next(); mr.Next(); }
  84.    tDelete(); gmx->ReleaseAndDelete(); return gmx;
  85. }
  86.  
  87. GeneralMatrix* ConstMatrix::Evaluate(MatrixType mt)
  88. {
  89.    if ((int)mt==0 || mt==cgm->Type()) { REPORT return (GeneralMatrix*)cgm; }
  90.  
  91.    REPORT
  92.    GeneralMatrix* gmx = cgm->Type().New(cgm->Nrows(),cgm->Ncols());
  93.    MatrixRow mr((GeneralMatrix*)cgm, LoadOnEntry);//assume won't change this
  94.    MatrixRow mrx(gmx, StoreOnExit+DirectPart);
  95.    int i=cgm->Nrows();
  96.    while (i--) { mrx.Copy(mr); mrx.Next(); mr.Next(); }
  97.    gmx->ReleaseAndDelete(); return gmx; // no tDelete
  98. }
  99.  
  100. GeneralMatrix* ShiftedMatrix::Evaluate(MatrixType mt)
  101. {
  102.    GeneralMatrix* gm=bm->Evaluate();
  103.    if ((int)mt==0)
  104.       { MatrixType mteqel(MatrixType::EqEl); mt = gm->Type()+mteqel; }
  105.    int nr=gm->Nrows(); int nc=gm->Ncols();
  106.    if (!(mt==gm->Type()))
  107.    {
  108.       REPORT
  109.       GeneralMatrix* gmx = mt.New(nr,nc);
  110.       MatrixRow mr(gm, LoadOnEntry); 
  111.       MatrixRow mrx(gmx, StoreOnExit+DirectPart);
  112.       while (nr--) { mrx.Add(mr,f); mrx.Next(); mr.Next(); }
  113.       gmx->ReleaseAndDelete(); gm->tDelete(); return gmx;
  114.    }
  115.    else if (gm->reuse()) { REPORT gm->Add(f); return gm; }
  116.    else
  117.    {
  118.       REPORT GeneralMatrix* gmy = gm->Type().New(nr,nc);
  119.       gmy->ReleaseAndDelete(); gmy->Add(gm,f); return gmy;
  120.    }
  121. }
  122.  
  123. GeneralMatrix* ScaledMatrix::Evaluate(MatrixType mt)
  124. {
  125.    GeneralMatrix* gm=bm->Evaluate();
  126.    int nr=gm->Nrows(); int nc=gm->Ncols();
  127.    if ((int)mt==0) mt = gm->Type();
  128.    if (mt==gm->Type())
  129.    {
  130.       if (gm->reuse()) { REPORT gm->Multiply(f); return gm; }
  131.       else
  132.       {
  133.          REPORT GeneralMatrix* gmx = gm->Type().New(nr,nc);
  134.          gmx->ReleaseAndDelete(); gmx->Multiply(gm,f); return gmx;
  135.       }
  136.    }
  137.    else
  138.    {
  139.       REPORT
  140.       GeneralMatrix* gmx = mt.New(nr,nc);
  141.       MatrixRow mr(gm, LoadOnEntry); 
  142.       MatrixRow mrx(gmx, StoreOnExit+DirectPart);
  143.       while (nr--) { mrx.Multiply(mr,f); mrx.Next(); mr.Next(); }
  144.       gmx->ReleaseAndDelete(); gm->tDelete(); return gmx;
  145.    }
  146. }
  147.  
  148. GeneralMatrix* NegatedMatrix::Evaluate(MatrixType mt)
  149. {
  150.    GeneralMatrix* gm=bm->Evaluate();
  151.    if ((int)mt==0) mt = gm->Type();
  152.    int nr=gm->Nrows(); int nc=gm->Ncols();
  153.    if (mt==gm->Type())
  154.    {
  155.       if (gm->reuse()) { REPORT gm->Negate(); return gm; }
  156.       else
  157.       {
  158.          REPORT
  159.          GeneralMatrix* gmx = gm->Type().New(nr,nc); gmx->ReleaseAndDelete();
  160.          gmx->Negate(gm); return gmx;
  161.       }
  162.    }
  163.    else
  164.    {
  165.       REPORT
  166.       GeneralMatrix* gmx = mt.New(nr,nc);
  167.       MatrixRow mr(gm, LoadOnEntry); 
  168.       MatrixRow mrx(gmx, StoreOnExit+DirectPart);
  169.       while (nr--) { mrx.Negate(mr); mrx.Next(); mr.Next(); }
  170.       gmx->ReleaseAndDelete(); gm->tDelete(); return gmx;
  171.    }
  172. }   
  173.  
  174. GeneralMatrix* TransposedMatrix::Evaluate(MatrixType mt)
  175. {
  176.    REPORT
  177.    GeneralMatrix* gm=bm->Evaluate();
  178.    if ((int)mt==0) mt = gm->Type().t();
  179.    GeneralMatrix* gmx=gm->Transpose(mt);
  180.    return gmx;
  181. }
  182.    
  183. GeneralMatrix* RowedMatrix::Evaluate(MatrixType mt)
  184. {
  185.    GeneralMatrix* gm = bm->Evaluate();
  186.    GeneralMatrix* gmx = new RowVector; MatrixErrorNoSpace(gmx);
  187.    gmx->nrows = 1; gmx->ncols = gmx->storage = gm->storage;
  188.    return gm->BorrowStore(gmx,mt);
  189. }
  190.  
  191. GeneralMatrix* ColedMatrix::Evaluate(MatrixType mt)
  192. {
  193.    GeneralMatrix* gm = bm->Evaluate();
  194.    GeneralMatrix* gmx = new ColumnVector; MatrixErrorNoSpace(gmx);
  195.    gmx->ncols = 1; gmx->nrows = gmx->storage = gm->storage;
  196.    return gm->BorrowStore(gmx,mt);
  197. }
  198.  
  199. GeneralMatrix* DiagedMatrix::Evaluate(MatrixType mt)
  200. {
  201.    GeneralMatrix* gm = bm->Evaluate();
  202.    GeneralMatrix* gmx = new DiagonalMatrix; MatrixErrorNoSpace(gmx);
  203.    gmx->nrows = gmx->ncols = gmx->storage = gm->storage;
  204.    return gm->BorrowStore(gmx,mt);
  205. }
  206.  
  207. GeneralMatrix* MatedMatrix::Evaluate(MatrixType mt)
  208. {
  209.    GeneralMatrix* gm = bm->Evaluate();
  210.    GeneralMatrix* gmx = new Matrix; MatrixErrorNoSpace(gmx);
  211.    gmx->nrows = nr; gmx->ncols = nc; gmx->storage = gm->storage;
  212.    if (nr*nc != gmx->storage)
  213.       MatrixError("Incompatible dimensions in CopyToMatrix");
  214.    return gm->BorrowStore(gmx,mt);
  215. }
  216.  
  217. GeneralMatrix* GetSubMatrix::Evaluate(MatrixType mt)
  218. {
  219.    REPORT
  220.    GeneralMatrix* gm = bm->Evaluate();
  221.    if (row_skip+row_number > gm->Nrows() || col_skip+col_number > gm->Ncols())
  222.       MatrixError("SubMatrix dimension error");
  223.    if (int(mt)==0) mt = MatrixType::Rect;
  224.    GeneralMatrix* gmx = mt.New(row_number, col_number);
  225.    int i = row_number;
  226.    MatrixRow mr(gm, LoadOnEntry, row_skip); 
  227.    MatrixRow mrx(gmx, StoreOnExit+DirectPart);
  228.    MatrixRowCol sub;
  229.    while (i--)
  230.    {
  231.       mr.SubRowCol(sub, col_skip, col_number);   // put values in sub
  232.       mrx.Copy(sub); mrx.Next(); mr.Next();
  233.    }
  234.    gmx->ReleaseAndDelete(); gm->tDelete(); return gmx;
  235. }   
  236.  
  237. GeneralMatrix* ReturnMatrix::Evaluate(MatrixType mt)
  238. { return gm->Evaluate(mt); }
  239.  
  240.  
  241. void GeneralMatrix::Add(GeneralMatrix* gm1, real f)
  242. {
  243.    REPORT
  244.    register real* s1=gm1->store;
  245.    register real* s=store; register int i=storage;
  246.    while (i--) *s++ = *s1++ + f;
  247. }
  248.    
  249. void GeneralMatrix::Add(real f)
  250. {
  251.    REPORT
  252.    register real* s=store; register int i=storage; while (i--) *s++ += f;
  253. }
  254.    
  255. void GeneralMatrix::Negate(GeneralMatrix* gm1)
  256. {
  257.    // change sign of elements
  258.    REPORT
  259.    real* s1=gm1->store; real* s=store; register int i=storage;
  260.    while(i--) *s++ = -(*s1++);
  261. }
  262.    
  263. void GeneralMatrix::Negate()
  264. {
  265.    REPORT
  266.    real* s=store; register int i=storage; while(i--) { *s = -(*s); s++; } 
  267. }
  268.    
  269. void GeneralMatrix::Multiply(GeneralMatrix* gm1, real f)
  270. {
  271.    REPORT
  272.    register real* s1=gm1->store;
  273.    register real* s=store; register int i=storage;
  274.    while (i--) *s++ = *s1++ * f;
  275. }
  276.    
  277. void GeneralMatrix::Multiply(real f)
  278. {
  279.    REPORT
  280.    register real* s=store; register int i=storage; while (i--) *s++ *= f;
  281. }
  282.    
  283.